% fields E for rho < a

function [rho, phi, z] = func_E_inside_check(r, m, CONSTS, coef, p1)

    k0    = CONSTS.k0;
    a     = CONSTS.a;
    eta   = CONSTS.eta;
    eps   = CONSTS.eps;
    g     = CONSTS.g;
    
    if (isnan(p1))
        rho = 0;
        phi = 0;
        z   = 0;
        return;
    end

    
    B1 = coef(1);
    B2 = coef(2);

    q1     = func_q1_from_p(p1, CONSTS);
    q2     = func_q2_from_p(p1, CONSTS);
    n1     = func_n1_from_p(p1, CONSTS);
    n2     = func_n2_from_p(p1, CONSTS);
    alfa_1 = func_alfa_1(p1, CONSTS);
    alfa_2 = func_alfa_2(p1, CONSTS);

    Q1 = k0*a*q1;
    Q2 = k0*a*q2;
    if(m>=0)
        bessel_m_q1 = besselj(m, Q1*r);
        bessel_mp1_q1 = besselj(m+1, Q1*r);
        bessel_m_q2 = besselj(m, Q2*r);
        bessel_mp1_q2 = besselj(m+1, Q2*r);
    else
        m1=abs(m);
        bessel_m_q1 = (-1)^m1*besselj(m1, Q1.*r);
        bessel_mp1_q1 = ((-1)^(m1-1))*((2*m1./Q1).*besselj(m1,Q1.*r) - besselj(m1+1, Q1.*r));
        bessel_m_q2 = (-1)^m1*besselj(m1, Q2.*r);
        bessel_mp1_q2 = ((-1)^(m1-1))*((2*m1./Q2).*besselj(m1,Q2.*r) - besselj(m1+1, Q2.*r));
    end

    rho  = -(B1*(((n1*p1+g)/eps)*bessel_mp1_q1 + alfa_1*m*(bessel_m_q1./(Q1*r))) + ...
                      B2*(((n2*p1+g)/eps)*bessel_mp1_q2 + alfa_2*m*(bessel_m_q2./(Q2*r))));
    
    phi  = 1i*(B1*(bessel_mp1_q1 + alfa_1*m*(bessel_m_q1./(Q1*r))) + ...
                        B2*(bessel_mp1_q2 + alfa_2*m*(bessel_m_q2./(Q2*r))));
    
    z    = (1i/eta)*(B1*n1*q1*bessel_m_q1 + B2*n2*q2*bessel_m_q2);
    
end